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Stripe formation in high- Tc superconductors 
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The non-uniform ground state of the two-dimensional three-band Hubbard model for the oxide high-T^ 
superconductors is investigated using a variational Mont Carlo method. We examine the effect produced 
by holes doped into the antiferromagnetic (AF) background in the underdoped region. It is shown that 
the AF state with spin modulations and stripes is stabilized due to holes travelling in the CuO plane. 
The structures of modulated AF spins are dependent upon the parameters used in the model. The effect 
of boundary conditions is reduced for large systems. We show that there is a region where incommensu- 
rability is proportional to the hole density. Our results give a consistent description of stripes observed 
by the neutron scattering experiments based on the three-band model for the CuO plane. 



I. Introduction 

A mechanism of superconductivity of high- Tc cuprates 
is not still clarified after the intensive efforts over a 
decade. An origin of the anomalous metallic properties 
in the underdoped region has also been investigated by 
many physicists as a challenging problem. In order to 
solve the mysteries of high- Tc cuprates, it is important to 
examine the ground state of the two-dimensional Cu02 
planes which are usually contained in the crystal struc- 
tures of high-Tc oxide superconductors.^ A basic model 
for the Cu02 plane is the two-dimensional three-band 
Hubbard model with d and p orbitals, which is expected 
to contain essential features of high-Tc cuprates. ^'^ The 
undopcd oxide compounds exhibit a rich structure of an- 
tiferromagnetic (AF) correlations over a wide range of 
temperature described by the two-dimensional quantum 
antiferromagnetism.''"* It is also considered that a small 
number of holes introduced by doping are responsible for 
the disappearance of long-range AF ordering. ^"^^ Recent 
neutron-scattering experiments have suggested an exis- 
tence of incommensurate ground states with modulation 
vectors given by = (vr ± 2tt5, tt) and Qc = (±47r(5, 0) 
(or Qs = (tt, TT ± 2ti5) and Qc = (0, ±A-k5) where 5 de- 
notes the hole-doping ratio. We can expect that the 
incommensurate correlations are induced by holes mov- 
ing around in the Cu-0 plane in the underdoped region. 

The purpose of this paper is to investigate the effect 
of hole doping in the ground state of the three-band 
Hubbard model in the underdoped region using a vari- 
ational Monte Carlo method^^"^® which is a tool to con- 
trol the correlation from weakly to strongly correlated re- 
gions. It is shown that the AF long-range ordering disap- 
pears due to extra holes doped into the two-dimensional 
plane. With respect to the initial indications given by 



the neutron-scattering measurements, the possibility of 
incommensurate stripe states is examined concerning any 
dependences on the hole density ^, especially regarding 
the region near of 1/8 doping. Although the possible 
incommensurate states are sensitively dependent upon 
the boundary conditions in small systems, the effect of 
boundary conditions is reduced for larger systems. 

The paper is arranged as follows. In Section II the wave 
functions and the method for the three-band Hubbard 
model are described. In Section III the results are shown 
and the last Section summarizes the study. 



II. 2D three-band Hubbard model and wave 
functions 

The three-band Hubbard model has been investigated 
intensively with respect to superconductivity (SC) in 
cuprate high-Tc materials. However, a non-uniform 
AF ground state for the three-band model has not yet 
been examined as intensively.'^^ The three-band Hubbard 
model is written as^^'"^^'^'^ 
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X and y represent unit vectors in the x and y directions, 

respectively, <t ^^"^ Pi±x/2.a denote the operators 

for the p electrons at the site Ri ± x/2, and in a similar 

way pJ±g/2,o- ^"^^ Vi±yli,a are defined. C/(= Vd) denotes 
the strength of Coulomb interaction between the d elec- 
trons. For simplicity wc neglect the Coulomb interaction 
among -p electrons. Other notations are standard and 
energies are measured in t^^p units. The number of cells 
which consist of d, Px and Py orbitals is denoted as N. 



The wave functions are given by the normal state, spin 
density wave (SDW) and modulated SDW wave func- 
tions with the Gutzwiller projection. For the three-band 
Hubbard model the wave functions for normal and SDW 



where Wxk = 2tdpSin(fca;/2), Wyk = 2tdpSm{ky/2), wu 
Kfe + ^Ikf'^ and Eu = [(e, - + <Y''- For 

the commensurate SDW Pk<T is given by a linear com- 
bination of dka, Pxka, Pyka, dk+Qa, Pxk+Qa and Pyk+Qa 

for Q = (tt, tt). Pg is the Gutzwiller projection opera- 
tor for the Cu d site. We can easily generalize it to the 
incommensurate case by diagonalizing the Hartrec-Fock 
Hamiltonian. The wave function with a stripe can be 
taken to be Gutzwiller, i.e. 

i^stripe = PG'tp°stripe- (6) 

ip'^tripe is tiic Slater determinant made from solutions of 
the Hartree-Fock Hamiltonian given as 

Htrial = H% + Y,[^nii - a{-lf'+y'mi]dtdi„, 

(7) 

where if^^ is the non-interacting part of the Hamiltonian 
H with the variational parameter and e^. The Slater 
determinant is constructed from wave functions of iV"e/2 
lowest eigenstates after diagonalizing H^riai in A;-space 
for each spin, where N^, is the number of electrons. 5ndi 
and TTii arc expressed by modulation vectors Qs and 
reprsenting the spin and charge part, respectively. In this 



states are written as 

i^n = PG n 4J0)> (2) 

|fe|<feF,<7 

^sdw = Pg n '^LlO)' (3) 

|fe|<feF,0- 

where aka is the linear combination of dka, Pxka and Pyka 
constructed to express an operator for the lowest band of 
a non-interacting Hamiltonian in the hole picture. Pq is 
the Gutzwiller operator given by 

PG = X{{l-{l-9)ndi^ndii), (4) 

i 

for Udia = c^i^dic- For tpp — 0, aka is expressed in terms 
of a variational parameter ip — id as follows: 



(5) 
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paper Sudi and rrii are assumed to have the form^^'^^ 

Sudi = -Yl «/cosh((xi - xf)/^c), (8) 
j 

rui = Aincom n tanh( (a;^ - xf (9) 
j 

with parameters a, Aincom, and ^g, where x**'" denote 
the position of a stripe. 

A Monte Carlo algorithm developed using the 
auxiliary-field quantum Monte Carlo calculations is em- 
ployed to evaluate the expected values for the wave func- 
tions shown above. ^^'•^^ Using the discrete Hubbard- 
Stratonovich transformation, the Gutzwiller factor is 
written as 

exp(-a^ndj|ndj4) = 

i 

^2^^ ^^M'^^Y^Siindit - Udii) - ^Yl^ridi-! + 

{si} i i 

ridii)], (10) 

where a = log(l/gf) and cosh(2a) = e"/^. The Hubbard- 
Stratonovich auxiliary field Si takes the values of ±1. The 
norm {tpstripe\'<P stripe) is Written as 



{ipstripeli/^stripe) = const. ^ JJ det(?!>o ' exp(F'^(u, a))exp(F'^(s, a))(?io ), (11) 

{Ui}{si} T 

I 

where V'^{s, a) is a diagonal 3N x 37V matrix correspond- 
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ing to the potential 

i i 

¥"{3,0) is given by ¥"{3,0) = diag(2a(TSi — 
a/2, • • • , 2aasN — a/2, 0, • • •) where diag(a, • • •) denotes a 
diagonal matrix with its elements given by the arguments 
a, V'^(s,a) has non-zero elements only for the d- 
electron part. The elements of {4>o)ij (i = 1, • • • , 3N;j = 
1, - ■ ■ ,Ne/2) are given by linear combinations of plane 
waves: 

{<l>oh = ""Mi^i ■ '*^)^« ' (13) 

for d-electron part (i = 1, • • • , N) where wfj is the weight 
of d electrons for £-th wave vector and j-th lowest level 
from below obtained from the diagonalization of H^^^ai- 
The p-elcctron parts arc similarly defined. Thus 

1,---,NJ2), (14) 

i^oh = ^exp(iri-k^)t«|,. {i = 2N+1, ■ ■ ■,3N,j = 

1,---,NJ2), (15) 

where wfj and w^- denote the weight of px and Py elec- 
trons, respectively. Then we can apply the standard 
Monte Carlo sampling method to evaluate the expecta- 
tion values. ^^'■^^ In order to perform a search for op- 
timized values of the parameters included in the wave 
functions, we employ a correlated-measurements method 
to reduce the cpu time needed to find the most descen- 
dent direction in the parameter space. "^^ In one Monte 
Carlo step all the Hubbard-Stratonovich variables arc up- 
dated once following the Metropolis algorithm. We per- 
form several 5 x 10^ Monte Carlo steps to evaluate the 
expectation values for optimized parameters. 

III. Antiferromagnetism and Stripes in the 
underdoped region 

We show the energy gain IS.Eaf for the uniform SDW 
state in reference to the normal state for optimized pa- 
rameters g, ip — id and AF order parameter Aa_f in Fig.l. 
The energy is lowered considerably by the AF long-range 
ordering up to about 20of U rs 8 — 12. 

IS.Eaf decreases monotonically as tpp increases and in- 
creases as U increases. One should note that AE^f is 
larger than the energy gain for the d-wave pairing state 
in the low-doping region near the doping ratio i5 ~ 0.1 by 
two order of magnitude. The boundary of the AF state 
in the plane of U and the hole density is shown in Fig. 2 
where AF denotes antiferromagnetic region. The doped 
holes arc responsible for reducing AF correlations which 
leads to an order-disorder transition. 

Let us look at doped systems on the two-dimensional 
plane with respect to modulated spin structures. Re- 
cent neutron-scattering measurements have revealed in- 
commensurate structures suggesting stripes. The 




Hole density 

FIG. 1. Uniform SDW energy gain per site with reference to 
the normal-state energy as a function of the hole density 5. 
Data arc from 8x8, 10 x 10, 12 x 12 and 16 x 12 systems 
for Sp — td = 2. For solid symbols U = 4 (circles), U = 8 
(squares), ?7 = 12 (triangles) and U = 20 (diamonds) for 
tpp = 0.2. For open squares U = 8 and tpp — and for open 
squares with slash U = 8 and tpp = 0.4. The lines are a guide 
to eyes. The Monte Carlo statistical errors are smaller than 
the size of symbols. 




Hole density 

FIG. 2. AF region in the plane of U and the hole density for 
tpp = 0.2 and tp — en = 2. 

AF states with spin modulations in space have been 
studied for the one-band Hubbard model'^^'''^"^^ and t-J 
model^^"^^ where various stripe structures are proposed. 
Our purpose is to examine the possible stripe structures 
and their parameter dependence based on the realistic 
three-band Hubbard model. We can introduce a stripe 
in the uniform spin density state so that doped holes oc- 
cupy new levels close to the original Fermi energy keeping 
the energy loss of AF background to a minimum. 

In the actual calculations we set = 1 and = 1 
in eqs.(8) and (9) since the expected values are mostly 
independent of and ^s- We optimize a in eq.(8) in- 
stead of fixing it in order to lower the expected energy 
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FIG. 4. Spin density {-lY''^Sz{e) (a) and hole density (b) 
functions at (5 = 1/8 where Sz{i) = natf — naii- Solid sym- 
bols axe for the 16 x 16 square lattice and open symbols are 
for the 16 X 4 rectangular lattice. The boundary conditions 
are antiperiodic in a;-direction and periodic in j/-direction, re- 
spectively. 



FIG. 3. Energy per site in reference to the normal state as a 
function of tpp for 16 x 4 and 16 x 16 lattices at 5 = 1/8. 
Circles, triangles and squares denote the energy for 4-lat- 
tice stripes, 8-lattice stripes, and commensurate SDW, respec- 
tively, where n-lattice stripe is the incommensurate state with 
one stripe per n ladders. In (a) the boundary conditions are 
antiperiodic in x-direction and periodic in j/-direction, and 
in (b) they are periodic in a::-dircction and antiperiodic in 
y-direction. The Monte Carlo statistical errors are within the 
size of symbols. 



FIG. 5. Spin structure in the incommensurate stripe state at 
J = 1/8. The boundary conditions are the same as in Fig.4. 
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value further because any eigenfunction of Htrial can be 
a variational wave funetion. It is also possible to assume 
that 5ndi and oscillate according to the cosine curves 
cos{An5xi) and cos(27r(5a;i), respectively, where 5 is the 
doping ratio. Both methods give almost the same results 
within Monte Carlo statistical errors. Let us define n- 
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FIG. 6. 



a function of tpp for 



Energy per site in reference to the normal state as 
a 16 X 16 square lattice at 5 = 1/16. 
Circles, triangles and squares denote the energy for 4-lat- 
tice stripes, 8-lattice stripes, and commensurate SDW, respec- 
tively. For solid symbols the boundary conditions are antiperi- 
odic in a;-direction and periodic in y-direction, and for open 
triangles they are periodic in x-direction and antiperiodic in 
jy-dircction, respectively. Monte Carlo statistical errors are 
smaller than the size of symbols. 
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lattice stripe as an incommensurate state with one stripe 
per n ladders for which the incommensurate wave vector 
is given by Qs = (tt ± ir/n,!:) and Qc — (±27r/n, 0) for 
the spin and charge parts, respectively. The incommensu- 
rate state predicted by neutron experiments at (5 = 1/8 
is four-lattice stripe for which = (tt ± 7r/4,7r) and 
Qc = (±7r/2, 0). In Fig. 3 we show the energy for commen- 
surate and incommensurate SDW states on the 16 x 16 
lattice at the doping ratio (5=1/8, where the incommensu- 
rability is given by 7r/4(= 27:5) for four-lattice stripes and 
7r/8 for eight-lattice stripes, respectively. The four-lattice 
stripe is stable in the range of 0.2 < tpp < 0.4. In Fig. 3 
we have shown the energy for two types of boundary con- 
ditions, which indicates that the effect of boundary con- 
ditions is not crucial for the 16 x 16 system, whilst the 
boundary conditions change the ground state completely 
for small systems such as a 16 x 4 lattice. The spin- 
correlation function exhibits an incommensurate struc- 
ture as shown in Fig. 4 and the hole-density function os- 
cillates corresponding to a formation of stripes. The 
spin structures are illustrated in Fig. 5. The energy at 
S = 1/16 is shown in Fig. 6 where the four-lattice stripe 
state has a higher energy level than for eight-lattice stripe 
for all values of tpp. The energy gain of the incommen- 
surate state per site in reference to the uniform AF state 
denoted as AEc-m is shown in Fig. 7 for tpp = 0.2, 0.25 
and 0.3. 



The incommensurability Ag/(27r) for tpp = 0.3 is also 
shown in Fig. 8 by solid circles, which is proportional 
to the doping ratio and is consistent with the neutron- 
scattering experiments for incommensurability.^" This 
should be compared with the variational Monte Carlo 
evaluations for the one-band Hubbard model^^ where the 
stripe states with large intervals are shown to be stable. 
Inorder to explain the linear dependence of Aq/{27r) on 
the hole density, the effect of tpp should be taken into 
account. The energy gain due to a formation of stripes is 
approximately proportional to the number of stripes. The 
size dependence of AEc-in is presented in Fig. 9; we ob- 
serve a tendency that AEc-in increases as the system size 
A'^ increases. The energy gain in the bulk limit is given by 
0m2tdp « 3meV for tpp = 0.3 where tdp = 1.5eV.^^"^^ 
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FIG. 7. Energy difference between the commensurate and in- 
commensurate states at (5 = 1/16 (16 x 16 lattice), 5 = 1/12 
(24 x 12 lattice) and 5 = 1/8 (16 x 16 lattice). From the top 
tpp = 0.2, tpp = 0.25 and tpp = 0.3. The boundary conditions 
are periodic in a-direction and antiperiodic in j/-direction, re- 
spectively. 



We present typical energy scales obtained from vari- 
ational Monte Carlo calculations in terms of t^p in Ta- 
ble I. The energy scales for superconductivity are con- 
sistent with experimental suggestions and energy differ- 
ence AEc-in between commensurate and incommensu- 
rate states are greater than the SC condensation energy 
by one order of magnitude. The commensurate AF en- 
ergy gain in reference to the normal state (denoted as 
AEaf) is larger than AEaf by one order of magnitude 
in the low- doping region. 
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TABLE I. Typical energy scales obtained from variationaJ Monte Carlo calculations for i7 = 8, tpp ~ 0.3 and ep — = 2. Aaf 
and AEaf denote the magnitude of AF order parameter and the AF energy gain compared to the normal state at half-filling 
5 = 0, respectively. Asc and AEsc represent the optimized SC order parameter and SC energy gain at 5 ~ 0.2, respectively. 
The last column indicates experimental suggestions. 

doping ratio EiicTgy{tdp) Exp. 

0.01~0.015(=15~ 20meV) 10~ 20nieV^'^''^'^ 

AEsc ~ 0.0005(=0.75meV)29.33 0.17~0.26meV".58 

Aaf S = ~ 0.6(=900meV) 

AEaf 6 = ~ 0.06(=90meV) 

Aaf (5 = 1/8 ~ 0.4 

Aincom (5=1/8 ~ 0.6 

AEc-in 6 = 1/8 ~ 0.002(=3meV) 
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Hole density 

FIG. 8. Solid circles denote incommensurability Aq/{2n) for 
tpp = 0.3 where the incommensurability is proportional to 
the hole density. For large tpp values the incommensurabil- 
ity equals zero as shown by the open circles. The boundary 
conditions are the same as in Fig. 7. 



IV. Summary 

Wc have presented our evaluations for the two- 
dimensional three-band Hubbard model using the varia- 
tional Monte Carlo method. We have examined an effect 
produced by holes doped into the AF state in the low- 
doping region. The boundary of AF phase is dependent 
on U as shown in the phase diagram in Fig. 2. The inho- 
mogcneous states with stripes arc stabilized due to hole 
doping so that the energy loss of the AF background is 
kept to a minimum with the kinetic-energy gain of holes 
compared to uniform (commensurate) AF state. In large 
systems the effect of boundary conditions is reduced in 
our evaluations. The distance between stripes is depen- 
dent upon the transfer integral tpp between oxygen or- 
bitals in the three-band model. There is a region where 
incommensurability is proportional to the doping ratio 
6 when S is small and the energy gain due to a stripe 
formation is approximately proportional to the number 
of stripes. A linearity of the incommensurability is con- 
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FIG. 9. Energy difference between the commensurate and 
incommensurate states as a function of 1/N for tpp = 0.3 (cir- 
cles) and tpp = 0.4 (triangles). Solid symbols arc for rectangu- 
lar lattices (16x4, 24x6, •• •), and open symbols are for square 
lattices (16 x 16, • • •) The incommensurate state is assumed 
to be the 4-lattice stripe state. The boundary conditions are 
periodic in a;-direction and antiperiodic in y-direction, respec- 
tively. 

sistent with the neutron-scattering measurements.'*'^ It 
is expected that the inhomogeneity plays an important 
role in the underdoped region with respect to anomalous 
metallic properties in high-Tg superconductors. We have 
also shown the typical energy scales obtained from vari- 
ational Monte Carlo calculations. It has been already 
established that the condensation energy AEsc and the 
magnitude of order parameter for superconductivity are 
in reasonable agreement with the experimental results. 
The energy gain due to AF ordering is larger than AEsc 
by about two order of magnitude and the energy dif- 
ference between the commensurate and incommensurate 
states is larger than AEsc by one order. The order of 
AF energy gain in reference to the normal state approxi- 
mately agrees with that for the t-J model. Our evalua- 
tions seem to overestimate the antiferromagnetic energy 
because of the simplicity of the Gutzwiller wave func- 
tions, which may give a starting point for more sophis- 
ticated evaluations such as Green function Monte Carlo 
approaches. 
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